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Abstract. One of the targets of the recently launched Fermi Gamma-ray Space 
Telescope is a diffuse gamma-ray background from dark-matter annihilation or decay 
in the Galactic halo. N-body simulations and theoretical arguments suggest that the 
dark matter in the Galactic halo may be clumped into substructure, rather than 
smoothly distributed. Here we propose the gamma-ray-flux probability distribution 
function (PDF) as a probe of substructure in the Galactic halo. We calculate this 
PDF for a phenomenological model of halo substructure and determine the regions 
of the substructure parameter space in which the PDF may be distinguished from 
the PDF for a smooth distribution of dark matter. In principle, the PDF allows a 
statistical detection of substructure, even if individual halos cannot be detected. It 
may also allow detection of substructure on the smallest microhalo mass scales, ~ M^, 
for weakly-interacting massive particles (WIMPs). Furthermore, it may also provide a 
method to measure the substructure mass function. However, an analysis that assumes 
a typical halo substructure model and a conservative estimate of the diffuse background 
suggests that the substructure PDF may not be detectable in the lifespan of Fermi in 
the specific case that the WIMP is a neutralino. Nevertheless, for a large range of 
substructure, WIMP annihilation, and diffuse background models, PDF analysis may 
provide a clear signature of substructure. 



PACS numbers: 95.85.Pw, 95.35. -l-d, 98.35.Gi 



The Gamma- Ray- Flux PDF from Galactic Halo Substructure 



2 



1. Introduction 

It has long been a goal of astrophysics and cosmology to determine the distribution and 
nature of the dark matter that populates our Galactic halo. Only more recently have 
we begun to focus on the possibility to detect substructures in the Galactic halo [1]. In 
hierarchical structure formation, small gravitationally bound dark-matter systems form 
first and then merge to form progressively more massive systems. In each stage, some of 
the earlier generations of structure may remain intact after merging, and so the Milky 
Way halo may contain substructures over a wide array of masses. Scaling arguments 
suggest that substructures may continue all the way down to the smallest mass scales 
at which there is primordial power [2l [3] , although the precise details may be uncertain 
[4j. If weakly- interacting massive particles (WIMPs) [5] make up the dark matter, the 
cutoff mass should be in the range 10"^ — 10^ [6], and if axions [7] make up the dark 
matter, it may be as small as lO""*^^ [8]. 

If WIMPs make up the dark matter, there may be several avenues toward detecting 
them. With the launch of the Fermi Gamma-ray Space Telescope (formerly GLAST) [9], 
however, there is now particular attention being paid to detection of energetic gamma 
rays from dark-matter annihilation in the Galactic halo (see, e.g., Ref. [10] and references 
therein). While the diffuse flux from such annihilations have been considered for a 
long time [H], the possibility to detect substructure, through angular variations in the 
background, is more recent [121 [131 [H] . If is possible that individual substructures may 
be resolved p5j. Proper motions of the smallest microhalos have also been considered 

[mini. 

In this paper, we propose the one-point gamma-ray-flux probability distribution 
function (PDF) as a probe of halo substructure. If dark matter is smoothly distributed, 
then the variation in the number of diffuse-background photons from one pixel to 
another should arise only from Poisson fluctuations. If, however, there is substructure, 
there will be additional flux variations from pixel to pixel. This may provide another 
route — an alternative to the angular two-point correlation function [121 IH] — to detect 
substructure statistically, especially for the very smallest microhalo mass scales|il It 
may also allow measurement of the substructure mass function, under certain model 
assumptions outlined below. 

We illustrate with a phenomenological model for Galactic substructure in which a 
fraction / of the halo is made of dark-matter microhalos with a power-law mass function 
(with a lower mass cutoff Mmin) and a constant mass-to-gamma-ray-luminosity ratio 
T = MminZ-^^min- The uext Section introduces this model and discusses the constraints 
from the Energetic Gamma Ray Experiment Telescope (EGRET) [19] to the parameter 
space. In Section [3l we calculate the flux PDF for this model and discuss the translation 
to a discrete distribution of counts in each Fermi pixel. We provide in Section H] numerical 
results for the PDF for an illustrative model. Section [5] determines the regions of the 

J It has been similarly suggested [H] that background fluctuations may be used to learn about the 
traditional astrophysical sources contributing to the diffuse background. 
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parameter space in which the PDF of substructure can be distinguished from that 
of a smoothly distributed background. In Section [6] we summarize and comment on 
additional steps that must be taken to implement this probe. 



2. Substructure/ Annihilation Models and EGRET Constraints 
2.1. Halo Model and Microhalo Mass Function 

We assume that a fraction / of the dark matter in the Galactic halo is composed of 
objects with a power-law mass function duh/dMh oc M^", independent of Galactocentric 
radius r. We shall take a = 2 in this work when evaluating numerical results, but our 
approach will hold in general. The mass function obeys the relation 

r^Iyn!.^ Aon 

fp{r)= dM,Mh—^{r,Mh)^{M,)nn{r), (1) 

where p{r) is the density profile of the Milky Way halo, Mmin and M^ax are the masses 
of the smallest and largest subhalos, and in the last equality we define the mean mass 
(Mh) as well as spatial number density nh{r) of subhalos. From Eq. ([1]) and the assumed 
shape of the mass function, we obtain 

''"^''^ M^i,ln(M^,./M„i,)' ^ ' 

{Mh) =M^i„ln(M„ax/M^m), (4) 

where in Eq. (j3]), we assumed Mmin Mmax- We use the NEW ^20J profile, 

(r/r^)(l + r/rsY 

where ps = 5.4 x 10^'^ Mq pc"^ is the characteristic density, and Vg = 21.7 kpc is the 
scale radius. The density is set to zero beyond a cutoff radius r^ = 10 r^, which is 
approximately the virial radius (i.e., the concentration parameter is c = r^i^/rs ~ 10). 
This normalizes the virial mass of the Milky Way halo to be lO^^M©, and gives 
po = 7 X 10~^ Mq pc~^ as the local density at the solar radius (ro = 8.5 kpc). 

Following other studies [21], we normalize the mass function by using the results 
of simulations [22j to fix the fraction of mass contained in high-mass microhalos. 
Specifically, we choose / such that 10% of the total mass of the halo is contained 
in microhalos of mass 10^ — IO^^Mq. We then extrapolate the power-law mass function 
found by the simulations down to a cutoff mass Mram below the simulation resolution; 
taking Mmax = IO^^Mq, / then becomes a function of Mmin: 

/(Mmin) = 0.10 log(Mmax/Mmin)/ log(Mmax/lO'M0) . (6) 

We shall suppress the argument Mmin when referring to / below. 

For a halo model with Mmin = M®, we find that approximately 52% of the total 
halo mass is contained in roughly 4.8 x 10^^ microhalos in the specified mass range. The 
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number density of microhalos in the solar neighborhood is about 34 pc~^. In this paper, 
we will examine a class of halo models in which Mmin is a free parameter, and falls in 
the range 10~^ — 10'^ Mq^ predicted by WIMP kinetic decoupling studies. 



2.2. Microhalo Annihilation Models 

Let us assume that the microhalos have NFW density profiles. The integrated number 
luminosity Lh from WIMP annihilation in an microhalo with NFW profile parameters 
rg, c, and ps is given by 

^^^^''^ f dVp'^a{c)Kp,M,. (7) 



ri o 

ml jh 

Here, is the integrated number of photons per annihilating particle, (cxf ) is the 
thermally-averaged annihilation cross-section multiplied by the relative velocity, and 
is the mass of the WIMP. In the second equality, 

^^'^ = 3(ln(l + c)-c/(l + c)) 

is a numerical factor resulting from the volume integral (with a dependence on c), and 
we have defined 

= (9) 

mr, mr, J dE 

A A 

Here, dN^/dE is the photon spectrum per annihilating particle. For the Galactic halo, 
using the NFW profile parameters defined in the previous section, we find that 

Lmw = 1-2 X lO^KMlpc-^ = 5.1 x 10^^ fsT GeV^ cm'^ (10) 

We now also assume that the integrated gamma-ray number luminosity Lh of each 
microhalo is proportional to its mass Mh, with constant mass-to-light ratio T = Mh/Lh- 
Then, the luminosity function is dn^/dLh = T{dnh/dMh). Note that throughout this 
paper, the luminosity is the number (not energy) of photons emitted per unit time; 
similarly, we deal with number fluxes (fluences) and intensities. 

These assumptions are consistent with the results of simulations, which indeed 
roughly find that oc M/j. In particular, Ref. [22] finds that 



3 X 10-12 ^ (11) 



Lmw Jmw p'^^^mw \^q, 
in the range of their simulation, which resolves subhalos down to Mh ~ 4 x lO^M©. We 
shall assume this relation holds down to the microhalo masses under discussion in this 
paper. Note that Eq. [TT] essentially relates the microhalo NFW profile parameters r^, 
c, and Ps (which may be complicated functions of mass) to those of the Galactic halo, 
which were stated in the previous Subsection. 

Combining Eqs. ITlfTTl we can now parameterize the magnitude of the annihilation 
signal by the parameter K (or equivalently, T-^ = 3.6 x 10"^ K Mq pc~^ = 
0.14 K GeV cm"^ or Lmin = T-^^Mmin, which are both proportional to K). Given 
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that our halo model and microhalo mass function were parameterized by Mmin, we see 
that our overall model has two parameters. We now discuss a constraint on this model, 
arising from an intensity limit observed by EGRET. 

2.3. EGRET Gonstramts 

The gamma-ray intensity Ih{'4') (units of photons cm^^ sec^^ sr^^) from microhalos 
along a line of sight at an angular separation ip from the Galactic center can be estimated 
as 

hW = ^ dl dL,L,—^{ril,^P),L,) 

/ 



4^^y dlp{r{l,^)), (12) 

where / is the distance along the line of sight; i.e., = Tq + — 2ro/cos'?/'. Compare 
Eq. [12] with the intensity Ig{4') from annihilation in the smooth component of the 
Galactic halo, which contains a fraction 1 — / of the total halo mass: 

IgW = ^^\~ / dl (r(/, ^)) . (13) 



Ait 

Note that Ih and Ig depend differently on p{r{l,ilj)), causing them to vary differently 
with -0. 

Current upper bounds to the diffuse gamma-ray background from EGRET place 
an upper limit on Ih + Ig- However, because of the lower energy range of EGRET, 
these upper limits apply to energies in the range 0.1 GeV < £^ < 10 GeV. Fermi will be 
more sensitive to photons with energies above 10 GeV (due to larger volume and better 
angular resolution at higher energies). For any given annihilation model, we are thus 
interested in the signal of gamma rays above 10 GeV, but must also check to see that 
the constraint in the lower energy range is obeyed. We see that we must examine the 
energy dependence of Ih + Ig, and hence the annihilation photon spectrum dN^/dE, 
in order to properly apply these constraints. We shall consider the annihilation photon 
spectrum in two different scenarios. 

In the first scenario, we assume that the WIMP is a neutralino, resulting in a photon 
spectrum per annihilating particle fit by an analytic approximation given by Ref. [23] : 

dN^ _ 1 0.42e-®^' 

ItE ~ m^xi-s + 0.00014' ^ 

where x = E/m^. For the neutralino particle properties, we choose typical values used 
in the literature. We set (erf) = 3 x 10~^^ cm^ sec~^, which reproduces the observed 
dark matter density if the WIMP is a thermal relic. We also choose = 85 GeV; 
this choice maximizes K for photon energies above 10 GeV, and hence maximizes the 
annihilation signal. 

Along these lines, when discussing annihilation signals in this scenario, we shall 
redefine the parameter K in all of the relevant preceding equations by using Eq. [H] in 
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Eq. [HI and integrating only over the energy range of interest to Fermi {E > lOGeV). 
Definitions for T, L/j, Ih, etc. in this energy range follow. With these values, we find 
the annihilation parameter for the neutralino model 

Kn = 4.2 X 10^8 pc^ sec"^ M'^ = 9.9 x 10"^^ cm^ sec"^ GeV-^ (15) 

By choosing these properties, we fix the annihilation parameter K; our model then 
depends only on the single parameter M^ain- Hereafter, we shall refer to our overall 
model in this scenario as the "neutralino model." 

To constrain this model, we rule out values of Mmin that result in intensities exceed- 
ing upper limits on the diffuse gamma-ray background found by EGRET (see Ref. |24j). 
That is, for a given Mmin, we require dIh/dE + dlc/dE < dIohs/dE over the EGRET 
energy range. Ref. [21] found that the gamma-ray background as observed by EGRET 
is roughly isotropic (after masking out the Galactic plane and center), and is suitably 
parameterized by dhhs/dE ^ 2.7 x 10"^ (E/6.5 GeV)-^-! cm'^ sec-^ sr"^ GeV"\ 
For our choice of neutralino properties, the relative shapes of the background and anni- 
hilation spectra are such that if the constraint dIh/dE + dlc/dE < dIohs/dE holds at 
6.5 GeV, then it is also satisfied over the entire energy range; thus, it suffices to check 
the constraint at this energy. We find that the constraint is satisfied for all M^[^ in the 
range 10~^ — 10^ predicted by kinetic decoupling studies. 

We plot in Fig.[T]the angular dependence of the microhalo intensity Ih above 10 GeV 
for a model in this scenario with fiducial cutoff mass Mmin = M®J§| We also plot Jg, 
the angular dependence of the intensity above 10 GeV from dark-matter annihilation 
from a smooth component containing 1 — / of the total Galactic halo mass, in order to 
show that it varies more rapidly with ip than the angular dependence of the gamma-ray 
intensity from substructure. We also show the intensity dIohs/dE of the gamma-ray 
background as measured by EGRET, integrated above 10 GeV. 

In the second scenario, we assume that WIMPs annihilate into monoenergetic 
gamma-rays of energy = 10 GeV. We leave Ke as a free parameter. We approximate 
the upper limit from Ref. [25] to the gamma-ray line intensity, averaged over the 10° x 10° 
region around the Galactic center, by 2 x 10~^ (i?/GeV)~^/^ cm~^ sec~^ sr^^ over the 
energy range 0.1 GeV < -E < 10 GeV (see also Ref. [26j), and we then derive an upper 
limit, 

/T^i < 10'' fi,h,i sec-i ^ 9 x lO"'' fi,h,E GeV~' sec'^ (16) 

Here, T^^ = 0.14 Ke GeV cm ^ is the light-to-mass ratio under the assumption 
of mono-energetic annihilation. Also, fi^h,E = Ih,E/{Ih,E + Ig,e + Id,E) < 1 is the 
fraction of the total gamma-ray intensity at 10 GeV from the Galactic center arising 
from annihilation in microhalos, and depends on the residual intensity I^^e from any 
astrophysical backgrounds that may not have been subtracted in Ref. |25j . 

§ Note that if p{r) oc as r ^ 0, the intensity is formaUy infinite at — 0. However, the flux 
from any finite-size window about the Galactic center involves an integral over the intensity, and the 
divergence of at -0 = is such that the flux is always finite. 
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Figure 1. The intensity Ih{ijj) above 10 GeV from microhalos as a function of the 
angle i/j the line of sight makes with the Galactic center (solid) for the fiducial photon 
spectrum given by Eq. [14] (the "neutralino model" ) and a cutoff mass of Mmin = Af® . 
We also plot (dashed) the angular variation of the intensity /g(V') above 10 GeV 
from dark-matter annihilation of a smooth component that contains 1 — / of the total 
halo mass. Note that the variation with -ijj of the gamma-ray flux from substructure 
is not as dramatic as that from annihilation in a smooth component. The intensity 
dIohs/dE integrated above 10 GeV of the gamma-ray background as measured by 
EGRET (dotted) is also indicated. 



In the case I^^e is negligible (i.e., the line intensity limit is saturated, with 
the observed intensity arising entirely from annihilation in the smooth halo and 
substructure), then fi^h,E depends only on the halo model and is a function of Mmin. 
Using Eqs. [12] and [131 calculation shows that a good estimate is given by fi^h,E ~ 
0.086 (Mmin/M®)~°°^^. For a given Mmin, the intensity limit then provides an upper 
bound on Ke', using Eq. [16] gives 

Ke<2x 10-29 (^n.n/M^r°-'' ,^3 ^^^-1 (^7) 
J (-'"mill ) 

Thus, in this scenario we shall consider models parameterized by Mmin and Ke, 
constrained by this limit. We shall refer to models in this scenario as "line models." 
Comparing Eq. [T7]to Eq. [15], we note that these line models may have much larger fluxes 
than the neutralino model. 

In Sections |3]]ll we shall discuss the neutralino model, using the assumed form of 
the photon spectrum to predict the PDF for this fiducial model. In Section [SI we shall 
examine how observation of the PDF may place constraints on the J^^; — Mmin parameter 
space for line models. 
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3. Calculation of the PDF 

The Fermi angular resolution at energies above 10 GeV is roughly 0.1°; throughout 
this paper we shall assume square pixels of solid angle (0.1°)^. This implies that the 
background flux will be measured in ~ 4 x 10® beams on the sky. One can then make a 
histogram of the number of counts in each beam. Our goal here is to make predictions 
for the shape P{F) for the distribution of these fluxes, under the assumption that these 
photons come from dark-matter annihilation in a clumpy Galactic halo. 

Although P{F) will in general be a function of the line-of-sight direction ip, we shall 
suppress this dependence in much of the presentation, reinserting it later when required 
for numerical results. We also refer to all probability distribution functions as P{x); the 
particular function under discussion should be clear from the argument x. 

If the population of sources has a flux-density distribution Pi{F), then the 
probability P{F) to see a total flux F (integrated over all sources in the beam) in 
a given beam is [27] 

P(F) = ^-i{e^(^^^^(^»-i)}. (18) 

Here jF{x} is the Fourier transform of x and JF^^ its inverse, and the flux-density 
distribution Pi{F) is normalized to / dF Pi{F) = 1. The quantity 

^ 'beam J / j,/ ,/2 



M^) = ^^/ dl'l'^piril',^)), (19) 

V^^h) Jo 

is the mean number of sources in each beam of solid angle f^beam (in sr). We reproduce 
in the Appendix the derivation of Eq. (fTS!) originally provided by Ref. [27] (see also 

I)- 



3.1. Derivation of Pi{F) 

The flrst step is thus to flnd the flux-density distribution Pi{F) for individual sources 
in the beam. This depends on the luminosity function and on the spatial distribution 
of microhalos. The luminosity function is P{Lh) oc L^", where again we take a = 2. 
The probability for an individual microhalo to be at a distance / along a line of sight ip 
is P{l,ip) oc p{r{l,ip)). We take a maximum cutoff at IcijP), corresponding to a cutoff 
radius = r{lc{4'),ip). 

We then flnd Pi{F) is given by 



Pi(F,^) = j dldLhP{l,i^)P{Lh)5 (^F- 



oc 



I 

''min[lc{ip),l{Luiax,F)] 



/ dl l'p{r{l, ^)) {l^F)-^e {AnPF - L^i,) e (L,,ax - 47rZ2F) 
Jo 



ocF~" / dll^-^''p{r{l,tlj)), (20) 

Jl{L^in,F) 

where the step functions enforce the cutoffs in P{Lh), and l{Li,F) = (Lj/47rF)^/^. 
Note also the implicit cutoff in Pi{F) for F < Lmin/47rZ^. Eq. ( l20l) can be evaluated 
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numerically for a given value of the parameter Lmin = T~^M^[^. The result is presented 
in Fig. [2] for the neutralino model, with the fiducial cutoff mass Mmin = M^. 

Note that Eq. ( 120|) yields the familiar Pi{F) oc (conventionally written as 

N{> S) oc ^-3/2) for a homo geneous spatial distribution of sources with a general 
luminosity function, if the condition /c(V') ^ ^(-^max, is satisfied over the range of F 
of interest. Under this condition, Pi{F) will also asymptote to F"^/^ at large F for 
a non-pathological spatial distribution. However, if lc{ip) < l{L^a.^, F) for values of F 
within the range of interest, then there will be a break in Pi{F); Pi{F) will tend to 
at F for which the second condition holds, and will then tend to F^^/^ at higher F. 

For the problem under discussion, values of Lmax in the interesting regions of 
parameter space are such that Pi (F) is negligible in the F"^/^ regime. Thus, the essential 
"large-F" dependence of Pi{F) will be 

3.2. Calculation of the Gounts Distribution 

The function P{F) gives the probability to observe a flux F from annihilation in 
substructure. Unlike the function P(Fbg) giving the probability to observe a flux 
Fbg from smoothly-distributed background sources (such as annihilation in the smooth 




Figure 2. The flux-density distribution Pi{F,tp ~ 90°), normalized to unity, for the 
flux from an individual microhalo drawn from a population of microhalos with: (1) 
an NFW spatial distribution and a luminosity function oc i^^, for a fiducial value 
of the minimum cutoff luminosity L^in (solid); (2) an NFW spatial distribution and 
uniform luminosity Ly^in (dashed); and (3) a homogeneous spatial distribution and 
uniform luminosity Lmin (dotted). Note that the first distribution follows a power- 
law of in the probabilistically observable range of F (following the power-law 
of the mass function), and then tends to at extremely large F. The last two 

distributions also tend to F^^^^. We have assumed the neutralino model and a cutoff 

mass Afniin = Af®. 



The Gamma- Ray- Flux PDF from Galactic Halo Substructure 



10 



component of the Galactic halo or other diffuse backgrounds), P{F) will not be a Poisson 
distribution. 

We will observe these fluxes in terms of photons beam~^ year^^, or similar units. 
However, the limits from EGRET require that the mean photon count per beam per 
year be less than one; thus, even those beams with the highest photon counts will only 
observe some small integer number of photons per year. It follows that we will need to 
discretize the continuous variables F and -Frg- Furthermore, emission of photons is a 
Poisson process. Thus, let the total number of photons measured in a given beam over 
an observation period T be C ^ E{F -\- F^q) G N; here E is the exposure in a beam 
given in units of cm^ sec, and is given hj E ^ AT, where A ^ 2000 cm^ is the area of 
the detector. 

The discrete probability distribution P{C) is then given by the sum of Poisson 
distributions with mean E{F + Fbg) weighted by P{F): 



Jo 

The shape of the discrete distribution P{C) is generally very similar to that of the 
continuous distribution P{EF) and is only slightly modified at the low end. 

4. Numerical Results 

Fig. [3] shows the results of numerical tabulation of the PDF P{F) for the neutralino 
model with a fiducial cutoff mass of Mmin = ^e- The PDF has a peak at low F and a 
power-law tail at high F. Here, only flux from dark matter annihilation in substructure 
and the smooth halo are included. Fig. H] adds an additional diffuse background with 
intensity Id = 10"'' cm~^ sec~^ sr~^ above 10 GeV; in this case, flux from annihilation 
only comprises a small fraction of the total observed flux, and the substructure PDF 
may be difficult for Fermi to detect within the mission lifetime. 

If the mass-to-light ratio T is increased (equivalently, if K is decreased) , with Mmin 
held fixed, then the photon flux decreases. The entire distribution is then scaled down 
along the F-axis. If M^[^ is reduced, with T held fixed, then the relative width of 
the peak of the PDF decreases. This behavior can be understood by considering the 
limit Mmin 0; in this case, we should expect to recover a smooth spatial distribution, 
resulting in a delta-function P{F). Note that this dependence on Mmin implies that the 
peak of the PDF must be resolved in order to measure Mmin; if the peak is obscured 
by an extraneous diffuse background, as in Fig. HI then Mmin may be a degenerate 
parameter. 

The distribution P{C) for discretized counts C is also plotted in Figs. [3] and HJ 
for a ten-year Fermi exposure. Also plotted is the Poisson distribution for a smoothly- 
distributed diffuse background of the same mean flux. As the Figures indicate, the 
large-F power-law tail of the PDF is qualitatively different than the exponential falloff 
of the Poisson distribution with F. Thus, detection of substructure amounts to detection 
of such a power-law tail. 




(21) 
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Figure 3. The angular-averaged discrete probability distribution function P{C) 
for the total photon number C in a given beam (circles), for the ncutralino model 
with a cutoff mass of Mmin = M^. Angular bins with widths of Aijj = 20° 
were used in the averaging. Only counts from annihilation in substructure and the 
smooth halo component have been included. An observation period of 10 years has 
been assumed. The angular-averaged continuous P{F) (solid) and a fitted Poisson 
distribution (squares) are also plotted for comparison. We have normalized to 
A^beam = 47r/(0.1°)^ w 4 X 10^, the number of beams at the angular resolution limit. 
Poisson error bars are also shown. 




Figure 4. The same as in Fig. [31 but with an arbitrary diffuse background with 
intensity Id — 10^^ cm^^ sec^^ sr^^ above 10 GeV added. This additional background 
adds a large Poisson-like feature to P{C) at low C, which obscures the substructure 
power-law tail. This suggests that the neutralino model may be just outside the range 
of P(D) analysis, if the diffuse background is indeed this large. 
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Figure 5. The same as in Fig. 21 but for a line model witli Ke = 10 Kn, an observation 
period of five years, and a background intensity Id,E = 10^^ cm~^ sec~^ sr^^ at 
10 GeV. Note that although the intensity from the diffuse background is still 
many orders of magnitude above the mean expected intensity from annihilation 
in substructure, the substructure power-law tail is detectable with high statistical 
significance. Furthermore, for these model parameters, there will be no detectable 
individual microhalos. Thus, P(D) analysis may be useful in detecting substructure 
even when individual point sources are not detected. 



Furthermore, the power-law tail of P{F) follows the power-law tail of Pi{F). This 
is simply because single bright sources dominate beams with high F. However, as 
discussed in Subsection 13.11 the power- law tail of Pi (F) in turn follows the power- law 
of the mass function. For example, the power-law tail in Fig. [3] indeed follows an 
dependence. Thus, P{F) not only provides a method of substructure detection; it can 
also reveal the substructure mass function. 

5. Detectability 

Figs. [3] and m are plotted for the neutralino model in which a fiducial value of K is chosen. 
This model predicts a mean flux far below the EGRET continuum limit; even Fermi may 
have to observe for a period of at least ten years in order to detect significant numbers 
of photons in beams in the substructure power-law tail. However, the constraint on 
line models given by Eq. [16] allows for choices of the parameters Ke and Mmin that 
result in mean fluxes much closer to the EGRET line intensity limit. Line models that 
saturate the limit will produce signals that could be easily detected by Fermi within a 
year. Fig. [5] plots the PDF for a fiducial line model. 

Of course, Fermi will also be sensitive to a range of line models predicting fluxes 
below the EGRET bound. However, for line models with mean fluxes below a certain 
level, the amplitude of the substructure power-law tail will so reduced that it will be 
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impossible to detect, as it was for the neutralino model in Fig. m In this Section, 
we determine the regions of the line model parameter space in which the PDF can be 
distinguished from the Poisson distribution expected for a completely smooth or diffuse 
background of the same mean flux, over an observation period of five years. Combined 
with the EGRET limit, this analysis will show the region of allowed parameter space 
that can be probed by study of the PDF. 

We determine the signal-to-noise with which the PDF P{C) can be distinguished 
from the Poisson distribution p{{C),C) with the same mean count rate (C). The null 
hypothesis of no substructure can be eliminated at the 3a level if S/N > 3, where 



Here, we label the angular bins by the central value of the bin tpi. The quantity 
A^beam,bin(^j) IS the number of beams contained in each bin, Cmax(V'«) is the highest 
count observed in each bin, and (C)^^ is the mean of the best-fit Poisson distribution 
in each bin. Eq. (122!) then quantifies the difference between the discrete probability 
distributions P{C,ip) and p((C)^,C), comparing the substructure PDF with the 
Poissonian distribution expected from a diffuse background (which may have angular 
dependence). In Fig. |6], we plot the regions of the Ke-M^^-^ parameter space in which 
the value of S/N indicates that substructure can be detected. Also plotted are the 
regions of the parameter space ruled out already by the current EGRET upper limit to 
the diffuse background. 

Note that in Fig. [6l models fall in the Poissonian regime when the substructure 
power-law tail at high C is obscured by the Poisson- like feature at low C; in this 
regime, P(D) analysis cannot be used. Since the Poisson- like feature arises from 
the flux from extraneous diffuse backgrounds (from sources other than dark matter 
annihilation in substructure and the smooth halo), the demarcation of the Poissonian 
regime is ultimately determined by the level of these backgrounds. In determining 
the Poissonian regime, we have conservatively assumed an arbitrary diffuse background 
of Id,E = 10^^ cm~^ sec~^ sr~^; in this case, the EGRET intensity limit is severely 
unsaturated, and only a small fraction of the observed diffuse background arises from 
dark matter annihilation in substructure and the smooth halo. In practice, the actual 
level of these diffuse backgrounds will determine the Poissonian regime, which may then 
cover a smaller region of parameter space than the conservative estimate presented in 



Of course, detection of a nontrivial PDF is also intimately related to the criteria 
for detection of point sources. The number of sources observed with flux greater than 




(22) 



and 




(23) 
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Figure 6. The KE-Mynin parameter space for the hne models. On the vertical axis, 
Ke is scaled by K^^, the annihilation parameter for the neutralino model given by 
Eq. [13 We indicate the region that is already ruled out by the EGRET upper limit 
to the diffuse background; where there will be > 1, > 10, and > 100 detectable point 
sources with flux greater than the five-year point-source sensitivity of Fermi; and where 
measurements of the flux PDF cannot be distinguished from a Poisson distribution, 
for an observation period of five years. Angular bins with widths of Aip — 20° were 
used in the calculation of S/N, and regions near the Galactic center (-0 < 30°) were 
masked. 



F is given by 

N{>F) =27, J d^sin^ j dU J dlfj^^ir{l,^/j),L^). (24) 

Examining this equation shows that A^(> F) is only weakly dependent on the cutoff mass 
Mmin. Furthermore, since the observed microhalos essentially comprise a volume-limited 
sample, for Ke in the range of interest A^(> F) oc Ke (at lower Ke, N{> F) oc K^'^ 
as expected for a flux-limited sample). Numerical calculation of dN{> F)/dM shows 
that the observed microhalos will predominantly be those of higher mass > 10^ Mq] 
although lower-mass microhalos are far more numerous, Fermi will not be sensitive 
enough to detect them individually p,7J. 

In certain regions of the parameter space for which Fig. [6]indicates a nontrivial PDF, 
substructure will be detectable via detection of individual microhalos, even without a 
detailed analysis of the PDF. We plot these regions, taking the Fermi five-year 5a point- 
source sensitivity at 10 GeV of F 2 x 10~^° cm~^ sec~^ (note that this sensitivity 
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assumes the same background level as in our determination of the Poissonian regime). 
The advantage of the full PDF, however, is that substructure can be detected even in 
regions of parameter space where individual microhalos elude detection. Measurement 
of the detailed shape of the PDF can also provide more information on the microhalo 
mass function and/or spatial distribution in the halo than would be obtained simply by 
point-source counts; e.g., the slope of the power-law tail in the PDF depends upon the 
slope of the mass function. 

Note that we could have done a similar analysis for a more general WIMP and 
substructure model in which some of the parameters (e.g., annihilation cross-section 
and spectrum, or subhalo concentration/boost parameters) were allowed to vary. The 
EGRET continuum constraint allows for a large range of such models. However, for 
simplicity we have only considered the parameter space of line models and the specific 
substructure model assumed in Section [21 

6. Conclusions and Comments 

We have proposed that the distribution of fluxes measured in individual Fermi pixels can 
be used to probe the existence of substructure in the Galactic halo to very small mass 
scales. By characterizing fluctuations in the diffuse gamma-ray background in this way, 
the existence of Galactic substructure may be inferred statistically even if individual 
halos cannot be detected. This statistical approach should be viewed as complementary 
to the use of an angular correlation function [121 113]. Since the PDF is a convolution of 
the microhalo mass function and spatial distribution, constraints to the parameters of 
these distributions may be obtained by measuring the PDF. 

The full PDF we have calculated may be useful even in situations where individual 
microhalos can be detected. For example, the flux in a pixel with a 3a excess which is 
interpreted as detection of a single point source may actually be due to several point 
sources; the probability that this is so may be inferred from the PDF. 

We have illustrated the PDF that results in a phenomenological model for 
substructure parameterized a microhalo mass cutoff Mmin, and a mass-to-light ratio 
T. This is almost certainly an oversimplification. In more realistic models, the mass 
function may differ from the particular power law we have assumed. The mass-to-light 
ratio may depend on the microhalo mass, and there may even be a spread of luminosities 
for each mass. The spatial distribution of microhalos may not trace the Galactic halo. 
Similarly, contributions to the PDF from astrophysical backgrounds (e.g., from cosmic- 
ray spallation or extragalactic sources) may need to be considered before a complete 
comparison of our model predictions with data can be made [29] . 

In our P(D) analysis, we did not consider the dependence of the angular resolution 
on the photon energy. Furthermore, we have also assumed here that each microhalo 
will fall within a single resolution element of Fermi. Taking into account the finite 
angular size of each microhalo will reduce the length of the power-law tails in the 
PDF, and will decrease the region of parameter space in which the PDF can probe 
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substructure. However, note that individual extended sources will also be more difficult 
to detect than point sources. A generalization of Eq. [2l] will give a smaller number of 
detectable extended sources; the corresponding hues in Fig. [6] will also shift upwards. 
Thus, there will still be an appreciable region of parameter space in which the PDF 
can be used to detect substructure even if individual sources cannot be detected. 
Moreover, a conservative rough estimate of the size of these microhalos can be found by 
approximating the microhalo mass density ph, assuming a formation redshift of z ~ 100 
and a concentration parameter of c ~ 1 A simple calculation then gives the angular 
size of the closest and most extended microhalos as 6* ~ ifPo/Ph)^^^ ~ 4°/^/^. Thus, if 
the beam size is increased such that the majority of extended microhalos fall within a 
single beam, then the point source P(D) formalism presented here is roughly valid. A 
more careful generalization may be required for comparison to data. 

We leave the inclusion of these additional levels of complication to future 
work. In addition to these future directions, one may also consider going further 
by combining the angular-correlation and PDF approaches. For example, the full 
two-point flux probability distribution function can be calculated and may provide 
additional observables with which to constrain the models or to distinguish a dark- 
matter background from other astrophysical backgrounds. Again, this is left for future 
investigation. 

Acknowledgments 

This work was supported by the Sherman Fairchild Foundation (SA), DoE DE-FG03- 
92-ER40701, and the Gordon and Betty Moore Foundation. 

Appendix: Derivation of P{F) 

Here we derive the relation between the flux-density distribution Pi{F) and the flux 
PDF P{F). Such a calculation is termed ^^P{D) analysis" in the literature, as it was 
flrst performed for observations of faint radio sources that produced "deflections" of the 
measuring apparatus. This P{D) analysis is useful in determining if an observed diffuse 
background is actually composed of numerous faint point sources. If this is the case, 
then there will be fluctuations in the diffuse signal from the random Poisson clustering 
of point sources in each beam. The shape of P{F) thus depends not only on Pi{F), but 
also the mean number /i of sources [Eq. (fT9l) ] in each beam. 

We wish to flnd the probability distribution for a total flux F in a beam, given that 
it is the sum F = Fi of the fluxes Fi from individual microhalos. Each of the Fi is 
a random variable with probability distribution Pi(Fj). Furthermore, the number k of 
fluxes Fi entering into the sum is itself a random variable given by a Poisson distribution 
with mean /i. Let us call Pk{F) the probability that k random variables Fi sum to F; 
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i.e., the probability that k microhalos emit a total flux F . Then 

oo 

P{F) = Y,pMP,{F), (25) 

k=0 

where p (/i, k) is a Poisson probability distribution for k with mean /x. 

It now remains to determine Pk{F). For /c = 0, it is clear that Po{F) = S{F); Pi{F) 
is given. For A; > 1, Pk{F) is given by 

/»oo poo I ^ \ ^ 

Pfc(F) = / rfFi . . . / rfPj n Pi (Fi) U(P - V p.). (26) 

The easiest way to compute Eq. fl26|l is to note that the Dirac delta function 
transforms the integral into a convolution [271 HH]- To see this, let us examine the 
integral for k = 2: 

poo poo 

P2{F)= dF, rfP2Pl(Pl)Pl(P2)(5(P-(Pi+P2)) 

Jo Jo 

POO 

= J dF,P,iF,)P,iF-F,) 

= (Pi*Pi)(P). (27) 

It follows that Pfc(P) = {Pk-i * Pi)(P); then by induction, Pfc(P) is given by Pi(P) 
convolved (or autocorrelated) with itself k times. Using the convolution theorem, it 
then follows that 



Pfc(P) =^-^{j-{Pi(P)}'} 



where JF denotes a Fourier transform. Note that Eq. fl28l) also holds for k = (and 
trivially for k = 1). 

Inserting Eq. (128|) into Eq. (!25l) and using the linearity of the inverse Fourier 
transform, we find 

oo _^ ^ 



k=0 



ii^HPiim' 



k\ 

. k=0 

= jr-^[e^iHPiiF)}-i)y (29) 

Eq. fl2^ gives the desired relation for P(P) in terms of Pi(P) and fi. Although the 
presence of the inverse Fourier transform prevents further analytic simplification in 
general, this expression can be computed numerically using fast Fourier transforms on 
a discretized Pi(P). 
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